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1. Introduction 

Lattice gases are interacting particle systems encountered in both equilibrium and non¬ 
equilibrium statistical mechanics. They are used as microscopic models for various 
physical and biological phenomena [1]. At large scales, considering macroscopic 
observables instead of the individual particles, these systems often evolve in time 
by deterministic hydrodynamic conservation laws. Understanding better fluctuations 
beyond the hydrodynamic behaviour is recognized as crucial in order to build a 
general theory for non-equilibrium phenomena ElElIl]. In many cases, the stochastic 
processes describing these fluctuations at large scale are independent of the details of 
the microscopic dynamics. This universal character of the fluctuations makes it very 
desirable to have exact expressions describing their statistics. This can be achieved by 
considering specihc microscopic models simple enough so that they may be solved. This 
approach was successfully used in the past for equilibrium statistical mechanics, the 
Ising model being a notable example. 

Another such model is the asymmetric simple exclusion process (ASEP) [51 El (TJ 
muM, whose dynamics breaks detailed balance and has thus a true non-equilibrium 
steady state at stationarity. ASEP is known to be integrable in the sense of quantum 
integrability, also called stochastic integrability HD] in the context of classical stochastic 
systems where convergence to a stationary state is ensured by the fact that the evolution 
operator is real valued, unlike in more traditional quantum integrable systems with 
unitary evolution where the issue of thermalization is still not completely settled. 

It is usually possible to diagonalize exactly the evolution operators of integrable 
models for hnite size systems using Bethe ansatz. For ASEP this leads, at least in 
principle, to exact expressions for the fluctuations. A technical problem is however 
to take the large scale limit of the hnite size, hnite time formulas, which is usually 
complicated as it involves delicate asymptotics of large determinants with entries written 
in terms of solutions of a large system of coupled polynomial equations of high degree. 
The situation simplihes enormously for the totally asymmetric simple exclusion process 
(TASEP), a special case of ASEP, for which some determinants can be computed 
explicitly, and the polynomial system of equations essentially decouples. 

We consider in this paper the one-dimensional TASEP on a ring of L sites. Each 
site is either empty or occupied by one classical particle. The dynamics consists of local 
hopping of the particles from one site i to the next i-|-l if the latter site is empty. Particles 
hop with rate 1, i.e. a particle has a probability dt to move in a small time interval 
dt. The dynamics conserves the total number of particles N, and the average density 
p = N/L is constant in time. A conhguration C of the system can be described by the 
occupation numbers of the sites rji G {0,1}, z = 1,..., L, where Pi = 1 means that site i 
is occupied and p* = 0 corresponds to an empty site. Equivalently, a conhguration can 
be specihed by the positions of the particles Xj, j = 1,..., N, 1 < xi < ... < xn < L. 

The state of the system can also be described by a height function Hi, i = 1,..., L 
in a mapping to an interface growth model. The mapping consists in evolving the initial 
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height “ Vi)y built from the initial occupation numbers rj^ of TASEP, by 

the following dynamics: each time a particle moves from site i to site i + 1 , Hi increases 
by 1. Extending the occupation numbers to a periodic function rji, i & Z of period L, 
the height Hi is also periodic of period L and verihes at all time Hi = Hi_i +]) — rji for 
any site i. 

We are interested in the (total, time-integrated) current Q, equal to the number 
of times a particle has moved anywhere in the system between time 0 and time T. 
This is a dynamical observable whose value can not be deduced from the knowledge 
of the positions of the particles in the system at time T only, but depends also on 
the history from an initial state Cq. It is however directly expressible from the height 
representation of TASEP as the difference between the hnal and the initial mean height. 

The generating function of the current {e^^)cQ^c, where the averaging is taken over 
all realizations of the process conditioned on starting from initial conhguration Cq at 
time 0 and ending in conhguration C at time T, obeys a (deformed) master equation |llj . 
In terms of the corresponding deformed Markov matrix M( 7 ), the generating function 
is equal to na 

,,0^ (C|e^«h)|Co) 

' {C|e™l»>|Co> ' ^ 

The denominator, called Z in the following, is the probability to observe the system in 
conhguration C at time T for an evolution starting in Cq at time 0. 

The problem is known to be integrable, as M{'y) closely resembles the quantum 
Hamiltonian of the XXZ spin chain. It allows an exact treatment using Bethe ansatz 
to diagonalize M{'y) and rewrite the generating function as a sum over normalized 
eigenstates 

. ( 2 ) 

r 

with Z ensuring that the generating function equals 1 at 7 = 0. For hnite systems, the 
eigenvalues and eigenvectors can be computed numerically very efficiently using exact 
Bethe ansatz formulas, which allows accurate evaluation of ([2]) or other observables such 
as average density prohle and current in a non-stationary setting [13], see also [HI for 
another approach based on an exact expression [H] for the propagator of periodic ASEP. 
Bethe ansatz also allows exact calculations in the thermodynamic limit of large L, N 
with hxed density p = A^/L, 0<p<l. This is especially true for periodic TASEP, 
for which the nearly decoupling structure of Bethe equations reduces enormously the 
complexity of the calculations. This has lead in the past to exact formulas for the 
spectral gap [m nzi [HI US] and large deviations of the current mm- 

In order to study the thermodynamic limit of ([2|), one needs to specify additionally 
how the hnal time and the initial and hnal conhgurations behave for large system 
size. The suitable scalings are known from KPZ universality [2T1 |22l |23l [2ll [25l l26] . 
whose name comes from the Kardar-Parisi-Zhang equation [27], and which describes 
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universal features of the statistics of fluctuations in various interface growth models, 
driven-diffusive systems and directed polymers in random media. KPZ universality is 
characterized by spatial correlations on the scale for large time. We consider the 
relaxation scale T ~ on which the correlation length saturates to the full system size 
L. Initial and hnal conditions are then chosen to be well described by smooth density 
prohles on the full range of the system. This is however not sufficient due to propagation 
of density ffuctuations around the system which hide the KPZ ffuctuations generated by 
the dynamics that we are interested in. Over times T S> T, these density ffuctuations 
move ballistically at the velocity 1 — 2 p. In order to correct for this, we take an initial 
conhguration Cq corresponding to a hxed density prohle pi and a hnal conhguration C 
described by a density prohle pf moving at velocity 1 — 2 p. The current huctuations are 
then dehned by subtracting from Q the deterministic part corresponding to the typical 
hydrodynamical evolution on the Euler time scale T ~ L of the macroscopic density 


prohle from Burgers’ equation, described in section |2] and Appendix A 


With the previously mentioned scalings for the various quantities, the large L limit 
of the summand of ([2]) can be performed explicitly for the special cases of unit step [2B] 
and hat initial and hnal conhgurations, giving exact formulas for the generating function 
and probability density of current huctuations. These exact results are extended to 
general step initial and hnal conhgurations with densities p+ and p_ by very accnrate 
extrapolation of high precision hnite size Bethe ansatz numerics. The main results are 
summarized in section [31 with some technical details about Bethe ansatz relegated to 
Appendix B 


Section 0] is hnally devoted to the special case of an evolution conditioned on hat 
initial and hnal states, for which the summation over eigenstates can be performed 
explicitly. It allows to extract the behaviour of current huctuations ft when the rescaled 
time t oc T /is small. With some proper dehnition of t fITTD and ft 11181) . one hnds 
the large deviations P{ft = ~ exp(—— H(m))) with some known constant 

C, and S dehned in flSip . This is the main resnlt of the paper. The rather technical 
saddle point analysis leading to it is carried out in [Appendix C 


2. Deterministic leading orders of the current and Burgers’ equation 

In this section, we snmmarize various known results about the deterministic evolution 
of the large scale density prohle of TASEP on times T ~ L from inviscid Burgers’ 
equation. We deduce from this the deterministic leading orders for the total current on 
times T ^ L. 


2.1. Hydrodynamic evolution: inviscid Burgers’ eguation 

From the stochastic microscopic dynamics of TASEP, the occupation number rji of site 
i evolves in time by 

^ = U) - 0i+.) (3) 
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with an instantaneous current ji = rji-fll — rii). At large scales, a deterministic evolution 
emerges at leading order for the density prohle p{x, r), obtained by averaging occupation 
numbers rji over sites i ~ xL. On the Euler time scale T = tL, the density profile evolves 
in time by a hyperbolic conservation law with one conserved quantity, the inviscid 
Burgers eguation 

drp + d^j = 0 , (4) 

with current-density relation 

J = P(1 - P) , (5) 

and initial condition p{x,0) = po{x) determined by the initial configuration of TASEP, 
see e.g. [2]. From time and space reversal in ([2]), Burgers’ equation also describes the 
macroscopic evolution for r < 0 of TASEP conditioned on ending at time T = 0 in a 
hnal conhguration corresponding to a density profile pi of average p: more precisely, the 
reversed prohle p{x, r) = p{l — x, —r) is the solution of Burgers’ equation with initial 
condition p{x, 0) = pi{l — x). 

The solution to Burgers’ equation (jl]) is only well dehned locally in time, even with 
smooth initial condition: after a hnite time, the solution p{x, r) develops shocks, i. e. 
discontinuities in x at some point 2 ; with a density lower on the left side of the shock 
X < z than on the right side x > z. Indeed, the characteristics x{t) such that p(a:(r),r) 
is constant in time {i.e. p{x{t),t) = po(xo) with xq = a:(0)) verify x'(r) = l—2po(xo). In 
an interval [xo,xi] where po decreases, it implies that the velocity of the characteristics 
starting at xq moves faster than the one starting at xi, which leads to the formation of a 
discontinuity. This makes dlj) ill-dehned since the motion of the shock can not be derived 
from Burgers’ equation. Unicity is recovered by imposing the additional constraint that 
the solution of (jl]) has to conserve the total density of particles p = fg dx p{x, r) since 
the number of particles is conserved in TASEP. This is equivalent to considering the 
viscosity solution of (jl]), obtained by taking the limit of vanishing viscosity z/ —)■ 0 in 
the solution of Burgers’ equation with the additional viscosity term nd^p in the right 
hand side. 


2.2. Integrated current and height function 


On the Euler time scale T = tL, the total current per site up to time T is equal at 
leading order in L to Q/L ~ LQ^jpo], with 

2r[po] = [ da [ dxj{x,a). (6) 


Jo Jo 

The instantaneous current j{x,T) is built from the current-density relation ([5]) with 
p{x,t) solution of (|1|) with initial condition po. Naively, the integral of j{x,T) with 
respect to x over the whole system is constant in time for the inviscid Burgers’ equation 
since drj{x, r) can be written as a derivative with respect to space as drj = {l—2p)drP = 
— (1 — 2p)dxj = — 2pfl/6. This argument breaks down after the formation of the 
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first shock since then the integration over space has to be done between shocks whose 
positions depend on time. 

As in the microscopic model, it is nseful to dehne a height function associated to 
the density prohle of the system by 

h{x,T) = ho{x) + / daj(x,a) , (7) 

Jo 

with initial height equal to 

ho{x)= [ dy{p-poiy)). (8) 

Jo 

This height function is equal to the large L limit of the microscopic height of the interface 
of TASEP Hi/L averaged over sites i ~ xL. Burgers’ integrated current (jH]) is then 
related to the height function by 

Qr = h(r) - ho , (9) 

with hnal and initial mean heights h(r) = dx h(x, r) and ho = dx ho(x). 

From dlj), the height function verihes d-rh = p(l — p) and dxh = p — p, which 
implies that h is solution of a deterministic KPZ equation without smoothing term 
drh = p(l - p) - {dxhfl - (1 - 2p)dxh. 

2.3. Large time evolution for smooth initial condition 

During the evolution, the number of shocks can increase when new shocks appear 
and decrease when consecutive shocks merge. With smooth initial density prohle, the 
number of shocks become constant at some value M > 1 at large time, generically 
M = 1, the density prohle converges to the hat prohle of density p, and the instantaneous 
current j(x, r) converges for large r to the stationary current 

J = p(l-p). (10) 

Burgers’ total current ([6]) is then approximatively equal to Qflpo] — Jt. We are 
interested in the corrections to this stationary value. They depend on the whole 
evolution of the density prohle between time 0 and time r, which involves in general the 
formation and merging of several shocks. 

At large times, the density prohle between consecutive shocks is approximatively 
given in the reference frame moving at the stationary speed of characteristics 1 — 2p 
(called the moving frame in the following) by ramps with negative slope of the form 

p(x +(1-2p)r,r) ~p-. (11) 

The position k, which corresponds to a density exactly equal to p, is located somewhere 
between the two shocks considered. 

In the generic case where only one shock remains at large enough time, its position 
is equal to z{t) ~ (1 — 2 p)r + k + | modulo 1 by conservation of the density. If M > 2 
shocks survive at large times and never merge afterwards, their positions in the moving 
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frame are equal to { k - + k +)/2 with and the positions at which the ramps on 
the left and on the right side of the shock have density p. 

The large time behaviour of the system is thus governed by the positions k at which 
the ramps have density p. Each number k, is the initial point of a characteristics of the 
partial differential equation (jl]) that never meets shocks and thus exists for all times. 
Such characteristics are called divides |29]. They have been dehned more generally for 
hyperbolic conservation laws with concave (or convex) current-density relation, of which 
Burgers’ equation (j4]) is the simplest non-trivial example. The initial points k of divides 
are the solutions of Po(k) = p such that dp (po(p) — p) < 0 for all x, see [21] theorem 
11.4.1. Equivalently, they are the locations of the global minima of the initial height 
prohle ho dehned in dHj): 

ho{n) = minho(x) . (12) 

X 

This is physically reasonable for TASEP since in the mapping to an interface growth 
model, the height only grows from local minima of the interface, and thus one expects 
that the large time behaviour is governed by the global minima of the initial interface. 

In the generic case where only one shock subsists at large times, k is unique and is 
equal to the position in the moving frame of the center of the ramp to which the density 
prohle converges. For non-generic initial prohles, ho can have M > 2 global minima, 
which leads at large time to the existence of M shocks. The special case of a hat initial 
density prohle po{x) = p corresponds to a situation with no shocks. 


2 . 4 . Burgers’ current at large time 


From the relation dxh = p —p, the height h{x, r) can be written as an integral over space 
with upper bound x. The constant of integration is obtained from h((l — 2 p)r -|- k, r) = 
jT+hfln), which follows from taking the derivative with respect to r of h((l— 2 p)r-|-K, r) 
and using the fact that characteristics starting from k have density p. One hnds the 
expansion 

h{x + (1 - 2p)r, r) ~ Jr -F ho{K) + , (13) 

4r 

for X inside the interval between two consecutive shocks corresponding to the ramp 
associated to a global minimum k, of Hq. Integrating with respect to x for each ramp, 
(1^ gives an expansion for Burgers’ current Qr- In the generic case where the global 
minimum of is unique, one hnds 

Qr[po] ^Jt + 7^[po] + ^ , (14) 


with 


7^[po] = minho(a;) — ho = — / dxx{po{K + x) — p) . 


(15) 


The quantity 7^[po] vanishes for a hat initial prohle po{x) = p. Furthermore, for any 
initial prohle, one has TZ[po\ < 0: the particles in TASEP move less easily on average 
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when the density prohle is not flat, which rednces the total integrated cnrrent. The 
expression ffT^ is checked in Appendix A for some simple piecewise linear initial density 
profiles by solving explicitly Bnrgers’ eqnation and calcnlating the cnrrent at finite r 
from (jH]). 

If ho has M >2 global minima, the term of order 1/r in ffT4|l is replaced by ^ 
with Afc the length of the interval for which the fc-th ramp has density larger than p. 
With flat initial condition po{x) = p, the cnrrent is exactly eqnal to Jr with no higher 
order correction. 


2.5. Deterministic current for TASEP conditioned on the initial and the final state 

Bnrgers’ eqnation describes the deterministic cnrrent for TASEP on a time scale T ^ L. 
On a longer time scale T ^ L, the macroscopic density profile stays flat for essentially 
all the evolntion, leading to a total cnrrent per site eqnal to Q/L = JT at leading 
order in L. Considering an evolntion conditioned to start at time 0 in a conflgnration 
corresponding to a fixed density profile pi and to end at time T in a conflgnration 
corresponding to a density profile pf in the moving frame, the first correction to the 
stationary valne of the cnrrent comes from time intervals with size of order L at the 
beginning and the end of the evolntion. It is expressed in terms of the qnantity TZ 
defined in ffT5|) as 

+ (K|pJ + K|p(])L . ( 16 ) 

with pflx) = pf(l — x). 

In the next section, we stndy the flnctnations of Q beyond the deterministic valne 
ffTbj) on the KPZ time scale T ~ nsing results from Bethe ansatz for specific initial 
and final states. 


3. Fluctuations 


On the KPZ time scale T ~ the density profile is typically equal to the constant 
profile p, except for small time intervals of duration ~ L at the beginning and at the end 
of the time range, where the density profile evolves from Burgers’ equation (jl]). From 
KPZ universality, height fluctuations in the moving frame have an amplitude x/I 

and are correlated on the spatial scale ~ L. We define the rescaled time 


yp(i -p) 


(17) 


3.1. Generating function 

We consider a fixed initial configuration Cq corresponding at large scale to the density 
profile Pi, and a final configuration C corresponding in the moving frame to the density 
profile Pf independent of T. All density profiles are periodic with periodicity 1. Based 






Current fluctuations and large deviations for TASEP on the relaxation scale 


9 


on the results of section [2] and on the scaling of height fluctuation in KPZ universality, 
we dehne current fluctuations as 

Q Qdet 


e* = 


(18) 


v^(r^L3/2 • 

with the deterministic value of the current Qdet given by ffTHj) . 

We are interested in the statistics of the random variable Q. We consider the 
generating function ([T]), ([2]) with fugacity 

s 


7 = 


v^(r^L3/2 ■ 

From KPZ universality, one expects that 


(19) 


( 20 ) 


has a hnite limit when L —)■ cx) with the scaling dm for T, and initial and hnal 
conhgurations corresponding to hxed density prohles in the reference frames described 
above. We dehne 


Gfls) = lim (e^^*)co^c • (21) 

L^oo 

The average over histories in the previous equation can be computed from the 
decomposition ([2]) over normalized eigenstates. One has 

p O'Qdet _ 

Gfls) = lim —^^e^®’'(^)(C|?^^( 7 ))(V'r( 7 )|Co) • (22) 

L^oo Z ^^ 

r 

with 

^ = 5^e^^’-(°)(C|V^,(0))(V^,(0)|Co) (23) 

r 

the probability to hnd the system in conhguration C at time T for an initial conhguration 
Co- Typical eigenvalues E of the Markov matrix M(0) scale as ~ L with Iie{E/L) < 0, 
see hgure [TJ Since the number of eigenvalues with a given value of EfL is of order 
exp(sL) j30] with hnite ’’entropy” s = 5{E/L), these typical eigenvalues have a vanishing 
contribution to Z. Extrapolating the small EfL behaviour s(E/L) r\-/ \ElL\fl^ [30] to 
eigenvalues E ~ a < 1 closer to the stationary eigenvalue 0, we observe that if 
Oi > —3/2, the contribution of the entropy exp(sL*^^’''^“^/®), s > 0 can not compensate 
the vanishingly small contribution of equal to exp(cL"’''^/^), Re e < 0. Furthermore, 
the eigenvalues with largest non-zero real part scale as [16]. Therefore, only the 

eigenstates whose eigenvalues have a real part scaling as contribute to (123|) . These 
eigenvalues correspond to the tip of the peak located at 0 in hgure [H The same kind 
of reasoning can presumably be used for non-zero 7 ~ too, for which the scalings 

for the entropy of eigenvalues should not be modihed. 

In the following, we call first eigenstates the inhnitely many eigenstates whose 
eigenvalue has a real part scaling as when setting s = 0 (for s 7 ^ 0 , the real 

part of the eigenvalues gains a term of order see dSID, but this term does not 

depend on the eigenstate and thus factors out of (1221) 1. From Bethe ansatz, each 
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Figure 1. Spectrum of TASEP at half-filling p = 1/2. The graph on the left represents 
in the complex plane the eigenvalues (divided by L) of the Markov matrix M(0) with 
N = 12 particies on L = 24 sites. Brighter coiors in the middie correspond to many 
eigenvaiues, darker coiors on the borders to fewer eigenvaiues. The graph on the 
right corresponds to the asymptotics (E71) Xr.( 27 rcr) at rescaied fugacity s = 0 of the 
eigenvaiues ciosest to 0. 


eigenstate is characterized by N pseudo-momenta kj, j = integers or half¬ 

integers depending on the parity of N. For the stationary state, the pseudo-momenta 
form a Fermi sea, k^ = j — {N A l)/2. The hrst eigenstates can be understood as 
particle/hole excitations over this Fermi sea [3T], corresponding to moving some pseudo¬ 
momenta with \k^\ < N/2 close to ±iV/2 to excited values with \kj\ > N/2, still close 
to ±A^/2. These excitations can be conveniently labelled by 4 hnite sets of positive 
half-integers and A^ representing respectively the positions of the hole and particle 
excitations on both sides of the Fermi sea, see hgure [2J Each creation of a hole on one 
side of the Fermi sea must be accompanied by the creation of a particle on the same 
side of the Fermi sea for the hrst eigenstates: any imbalance leads to eigenvalues with 
real part scaling as with some a > —3/2. It implies that the cardinals of the sets 
verify the constraints 

= |A([| = |A+| and = jAg | = \A~\ . (24) 

In the following, the four sets A^ are collectively denoted by the index r. 


3.2. Large L asymptotics 


For each hrst eigenstate r, it is convenient to introduce the function Xr, with branch 
cuts i[7r, oo) and — i[7r, oo), dehned by 


87r3/2 




3 1 iu 


3 1 iu 








Current fluctuations and large deviations for TASEP on the relaxation scale 


11 



































































































































5 3 

2 2 

3 7 

2 2 


9 7 1 

2 2 2 

'-V-' 

1 5 9 

2 2 2 

'-v-' 


A- Ao A+ A+ 


Figure 2. Graphical representation of the (half-)integers kj (coloured squares) 
characterizing some of the first eigenstates, whose eigenvalue has a real part scaling 
as L~^P. From top to bottom, the pictures correspond to the stationary state, the 
two eigenstates giving the spectral gap, and a generic eigenstate with particle-hole 
excitations described by four sets of half-integers Af , . 


87r3/2 / 


2 eA+ 


m \ 3/2 




-1 a -b 




aGA- 

iu \ 3/2 


aeA+ 


1 a -b 


(25) 


m \ 3/2 


where ( is the Hurwitz zeta function. For the stationary state, the four sets are 
empty, and the function reduces to a polylogarithm from Jonquiere’s identity: Xo(u) = 
-(27r)"^/2Lij.^2(-e“) if —TT < Imu < TT. We also introduce the complex number 
Cr = Cr(s), solution of 


= s . (26) 

Neither existence nor unicity of has been proved; numerics seem however to indicate 
that both hold for any choice of the sets satisfying fl 2 ^ if the rescaled fugacity s verihes 
Res > 0, which is a consequence of the restriction Re 7 > 0 of Appendix B The 
stationary state has the singular solution Co(s) -A —00 when s ^ 0. For some choices 
of the sets A^, the solution Cr can have very large imaginary part, which makes 
an analytic continuation needed if one wants to work with polylogarithms instead of 
Hurwitz C functions. 

The eigenvalues of the hrst eigenstates have the large L expansion 


Erix) - 


p) 2m{l - 2 p)jx ^ 


P) 


XriflTlCr) , 


(27) 


x/I L L3/2 

see hgure [U for a graphical representation of the hrst few Xrifl'^Cr)- The special case 
of the spectral gap, corresponding to the hrst non-zero eigenvalue, was obtained in 
[TBl [T71 [TKl [T^ . It has been also studied for periodic ASEP [32], for TASEP [331 [31] 
and ASEP [351136] on an open interval, and for periodic ASEP with several species of 
particles [33 ES]. 
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We consider unnormalized Bethe eigenvectors, described more precisely in 
Appendix B[ The left and right eigenvectors cfr can be chosen in such a way that 




(28) 


where the conhguration x with particle at positions 1 < Xi < ... < ^ L and the 

conhguration x with particle at positions 1 < Xi < ... < x^ < L are related by space 
reversal Xj = L + 1 — xn+i-j- The large L limit for the normalization of these Bethe 
eigenstates has been obtained in [2B] : 


n 

_ 

(0r |0r) 



(29) 


with n = the total number of conhgurations. We changed the overall normalization 
of the eigenstates from [2H] in order to make the elements of the eigenvectors simpler, 
Appendix B 


see 


For a conhguration C corresponding to a hxed density prohle po, we write the 
asymptotics of the elements of the eigenvectors as 


(MC) - . 


(30) 


Shifting a conhguration by a distance X gives an additional factor |;o (C|0r); 

with in particular X = (1 — 2p)T for a conhguration corresponding to a density prohle 
hxed in the moving frame. Gathering everything, this implies for the generating function 
of current huctnations 


Gfls) 


\ g27rc,. gtXr(27r Cr) 

^ \/^x;'(27rc^) 




(31) 


with normalization constant Zt = VtZ equal to the probability of having the system in 
conhguration C at time T starting in Co, divided by the stationary probability 12“^. The 
stationary eigenvector r = 0 at fugacity 7 = 0 verihes <ho[p] = 1 independently of p, and 
Zt ^ 1 when f —)■ oo since Co(s) —)■ —oo when s —)■ 0. 

From section |2l the quantity <Fr[po] is expected to be independent of L and to 
depend only on po and not on the details of the conhguration C. It can be computed 
explicitly in the special cases where the Bethe ansatz expressions for the eigenvectors 
reduce to Vandermonde determinants. This is in particular the case for hat and unit 
step conhgurations, which leads to exact formulas for the current huctuations in four 
cases, denoted hat —)■ hat, step —)■ hat, hat —)■ step, step —)■ step, depending on the 
initial and the hnal state on which the evolution is conditioned. The hat —)■ hat case is 
studied in much detail in section HJ 


3.2.1. Flat configurations The component of the eigenvector is computed by 
elementary manipulations in Appendix B for a hat conhguration 3F with particles at 
positions Xj = X + (j — l)/p, j = 1,..., iV and p“^ integer, which corresponds at large 
scale to a hat prohle p^{x) = p. One has 7?.[p'^] = 0 and 

4,17] = 4f 


(32) 
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independently of X, with 


= 1 


{yl+=A-}-^{A(, =A+} _)_ g27rc^^l/4 ’ 


(33) 


The constraint Af = A~, Aq = A~^ implies = mf + m~. The elements of the 
eigenvector corresponding to flat conhgnrations vanish exactly when Aq A~ or 
Aq 7 ^ A~^ becanse of the symmetries of the conhguration. 

The factor (1 + in fl33|) is nnderstood with the same branch cnts i[|, oo) 

and —i[|, cxo) as Xr(27rc). With the nsnal dehnition of the non-integer power = 
exp(l log 2 ;) and the nsnal branch cnt for the logarithm, the factor (1 -|- is 

interpreted as (—i)d™h+i/2)J _|_ g27rc^-i/4 ^i^p the largest integer lower than x. 

The expression fl33|) has been checked nnmerically using rational Richardson 
extrapolation [32] (also called the Bulirsch-Stoer method) of hnite size Bethe ansatz 
numerics. Richardson extrapolation allows to extract the constant term /o of an 
expansion of the form f{L) = YlT=o knowing a few values f{L) (with high 

precision) for moderate values of L. It often allows to extract around 10 correct digits 
of /(O) knowing 20 values f{L) with only one signihcant digit in common with /o, see 
table [T] for an example. Richardson extrapolation comes naturally with an accurate 
estimator for the error on /q. It was used here not only for the conhguration X, but 
also for more general conhgnrations corresponding to a macroscopic hat prohle, built 
by repeating clusters of the form l^’+O'’” with r+,r_ > 0, r+ -|- r_ < 5. A perfect 
agreement was found with fl33|l within at least 10 digits, see tabled] for an example. All 
the computations were done with a generic value s = 0.2 -|- i for the rescaled asymmetry. 
The natural exponent ca = 1/2 was used for the extrapolation. 


3.2.2. Step configurations In the case of unit step conhgnrations Sx/l-, where sites 
from X to X -|- A^ — 1 are occupied while the rest of the system is empty, the density 
prohle is called Px/l- corresponding element of the eigenvector is a Vandermonde 
determinant. The calculation of its large L asymptotics is signihcantly more involved 
[28] than in the hat case, and can be obtained using two-dimensional Euler-Maclaurin 
formula in a triangular domain with combinations of logarithmic and square root 
singularities at all the edges and corners. From flTS]) . one has 7?.[pf] = —p(l — p)/2, 
and the eigenvectors are (see Appendix B ) 

(34) 


Appendix B 


and *hr[pf] = for the reversed prohle pf according to fl25]) . The quantity is 

equal to 




(2i) 


■ u/{A^)u/{A^)u{A+)u{A~)u{A^, A^)u{A+, A' 


X exp ( lim -ml log A + / du 

\ A —j ^ 2 / 




(35) 
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L 

Numerical value 

Richardson extrapolation 

4 

0.646361 - 0.409949 i 

o 

1 

o 

8 

0.658283 - 0.403781 i 

0.7-0.4i 

12 

0.66433 - 0.398974 i 

0.7-0.3i 

16 

0.668107 - 0.39556 i 

0.7-0.4i 

20 

0.670745 - 0.393016 i 

0.69 - 0.37i 

24 

0.67272 - 0.391037 i 

0.694- 0.367 i 

28 

0.674269 - 0.389443 i 

0.6941 - 0.3667 i 

32 

0.675525 - 0.388125 i 

0.6940 - 0.3666 i 

36 

0.67657 - 0.387012 i 

0.6941 - 0.3667i 

40 

0.677457 - 0.386056 i 

0.69407-0.36669 i 

44 

0.678222 - 0.385223 i 

0.694065 - 0.366690 i 

48 

0.678891 - 0.384489 i 

0.6940651 - 0.3666903 i 

52 

0.679482 - 0.383836 i 

0.6940651 - 0.3666903 i 

56 

0.68001 - 0.383251 i 

0.6940651 - 0.3666903 i 

60 

0.680484 - 0.382721 i 

0.694065124 - 0.366690292i 

64 

0.680914 - 0.382239 i 

0.694065123 - 0.366690292i 

68 

0.681305 - 0.381799 i 

0.6940651235 - 0.3666902918i 

72 

0.681663 - 0.381394 i 

0.6940651235 - 0.3666902918 i 

76 

0.681993 - 0.38102 i 

0.69406512350 - 0.36669029178 i 

80 

0.682298 - 0.380673 i 

0.694065123499 - 0.366690291776 i 

84 

0.682582 - 0.38035 i 

0.6940651234993 - 0.3666902917760 i 

88 

0.682845 - 0.380049 i 

0.6940651234993 - 0.3666902917760i 

92 

0.683092 - 0.379766 i 

0.6940651234993 - 0.3666902917760i 

96 

0.683323 - 0.379501 i 

0.694065123499272 - 0.366690291775979i 

100 

0.68354 - 0.379252 i 

0.694065123499272 - 0.366690291775980i 


Table 1. Richardson extrapolation of finite size Bethe ansatz numerics for the 
component of the right stationary eigenvector corresponding to a configuration C of 

the form _where • corresponds to a particle and _ to an empty 

site. The calculations are done with 50 digit precision and the rescaled fugacity is 
equal to the generic value s = 0.2 + i. The first column corresponds to the size L of 
the system and the second column to (the first digits of) the numerical evaluation of 
{C\(j)o) from the exact Bethe ansatz expression (IB.71) with the rescaling above (IB.12|) . 
The third column is the result of Richardson extrapolation with exponent w = 1/2 
from the numerical values for system size < L, truncated at the error estimated by the 
extrapolation method. The exact asymptotics, given by (ESI) with four empty sets, is 
equal within 20 digits to 0.69406512349927191436 - 0.36669029177597961516 i. 
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with combinatorial factors 

u}{A) = (a — a') and oj{A, A') = (a + a') . (36) 

a,a'&A a&Aa'&A' 

a>a' 

More generally, Richardson extrapolation of finite size Bethe ansatz numerics indicate 
that for any configuration 5^"’^+ corresponding at large scale to a step density prohle 
po with densities between x and x + a and p_ elsewhere, p_ < p+, one has 

= Q-^i^PrAPo]^s ^ ^ 37 ^ 

where k[po] = x + a is equal to the position of the center of the ramp in the moving frame 
at large times, defined more generally by flT^ . This was checked for configurations of 
the form for all the cases with 0<r_<r+<r<5 and 

£+ = £_, and for some other cases with £+ = 2£_ and £+ = i-/2, r_ < r+. An exact 
match was found with fl37j) within the error estimator of the extrapolation method. 


3.2.3. Generic configurations From section |2l any generic smooth density prohle leads 
asymptotically to the same linear decreasing prohle fllip for large ty/L, and the current 
huctuations are expected to be the same as in the step case. Checking this with high 
precision using Richardson extrapolation does not seem possible, however, due to the 
lack of a natural sequence of conhgurations Cl leading to a clean expansion in powers 
of l/y/L for the eigenvectors. Nevertheless, limited numerics on linear and sinusoidal 
prohles seem to conhrm the asymptotics fl3n|l with TZ[po] given by flT5|) . These numerics 
also seem to indicate the presence of extra non-universal constants shifting the current, 
that have to be removed in order to recover fl37)) . The non-universal constants seem to 
vanish for the left eigenvector with linear increasing prohle and the right eigenvector 
with linear decreasing prohle, which is probably related to the fact that the expression 
flTT]l for the density prohle at large time is exact for a linear increasing initial density 
prohle, as in the case of a step prohle, see [Appendix A 


3.3. Probability distribution of the current fluctuations 


The probability density Pt of the current huctuations ft is obtained from Fourier 
transform of the generating function fIM]) . One has 

ds 

Ptiu)= — (38) 

J-oo 

The probability density is plotted in the hat —)■ hat case in hgure |3l based on numerical 
evaluations where only a hnite number of eigenstates r are kept in flTT]) . More eigenstates 
are needed to ensure reasonable convergence for small values of t. 

Making the change of variables s —?■ c removes the necessity to solve (|26|) in the 
expression fl55]) . (HU]) for the probability density. One hnds 


T/3r 


P,(n) = 


Z, 


t Je 1 


dc 


.,27rc 




(39) 
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Figure 3. Plots as a function of u of the probability density PP of rescaled 
current fluctuations computed numerically from (1551) . (HUl) in the flat -r flat 

case. The integral over s in (I38p is discretized as a sum over 101 equally spaced values 
between —5/t and 5/t. The generating function G is evaluated from (1401) by summing 
over the 12 eigenstates corresponding to sets A = Af = A~ and A = Aq = A~^ with 
Jfia&A ® + '^aGA ^ The three graphs correspond respectively to rescaled time 

t = 0.2 (flattest curve), t = 0.1 and t = 0.05 (most peaked curve). When t —>■ 0, the 
random variable converges with probability 1 to u* « 1.291468, indicated by a 

dotted line. 


with Cr = Cr(is) replaced by c in the expressions for The integration range follows 
from the large |s| asymptotics 2 TTCr{s) ~ ^ ^ The quantities Cr(is), 

s G M are plotted in figure 0] for some eigenstates. 

4. Large deviations in the fiat —flat case 

In this section, we specialize to the case of an evolution conditioned on flat initial and 
flnal conhgurations. There, the sum over eigenstates in the probability density of current 
fluctuations can be computed explicitly as an inhnite product. From this expression, an 
exact formula is derived for the large deviations of current fluctuations at short rescaled 
time t. 

4-1. Probability distribution 

If both pi and p/ correspond to flat profiles of density p, the generating function (IMll 
simplihes to 

\ g27rCrgtXr-(27rCr) 

1 {^0 -^+1 Q^-KCr ■^> 4 { 2 TTCr ) 


Gt{s) 


(40) 
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Figure 4. Plots in the complex plane of the curves c,.(is), s £ M for the few first 
eigenstates r contributing to the current fluctuations in the flat —>■ flat case. The curve 
for the stationary state r = 0 goes to —oo when s —>■ 0. The curves for all the other 
eigenstates verify Reci.(is) > 0. When s —>■ ±oo, the curves diverge to e^‘’^/^oo. The 
thick, red lines correspond to the branch cuts i[|, oo) and —i[i, oo) in the variable c of 
the integrand in (H^ . 


The probability density of current fluctuations fl3^ is then equal to 


g27rCg-tix[,(27rc)gtxo(27rc) 


P,(«) ^ ^ 


e 1 


\\P2t: \/l + 


i27rc 


E 1 




^ VIV^) 


(41) 


X e 


-^(E vd(a-ic)3/^+E v^(a+ic)3/^) 


where Xo is the function (I2S]) for the stationary state, corresponding to four empty sets. 
The sum over the sets can be computed by adding a contour integral to enforce the 
constraint |Aq | = |Ag |. One has 


ir/3^ 


Pt(u) = 


(dc g 27 rCg-'UXo( 27 rc)gtxo( 27 rc) r 




n [( 


\/2^ vTTe^ I 2i7rz 

3l^_i6,A^EE(a-ic)3/2 


/ —WTTUyJ a- 


aGN+J 


1 + ze 


(42) 


f.2. First cumulants of the current 

The hrst cumulants of ft are obtained by taking derivatives with respect to s at s = 0 of 
the generating function Gfls) (ITUD . They can be computed numerically by truncating 
the sum over hrst eigenstates to keep only the eigenstates corresponding to small values 
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Figure 5. Rescaled first cumulants of the current fluctuations plotted 

as a function of t, for an evolution conditioned on flat initial and final states. The solid 
lines correspond to values of the cumulants obtained by taking derivatives at s = 0 
of the generating function (HUl) . computed by summing over the 6639349 eigenstates 
corresponding to sets A = Af = A~ and A = A^ = A+ with ® + SaeA “ — 

From top to bottom, they represent the rescaled average (blue), variance (red), fourth 
cumulant (orange) and third cumulant (cyan). The dashed lines represent the short 
and long time values, given respectively by (IST)) and (1171) . 


of XlaeAj ® + XlaeAp ® + X]aeA+ ® + XlasA" ® eigenstates converges 

rather quickly if the rescaled time t is not too small. This is especially true when either 
the initial or the hnal state is flat: then most eigenstates do not contribute because 
of the constraints Aq = A~ and Aq = A'^ on the sets. The four hrst cumulants are 
plotted in the flat -A flat case as a function of time in hgure [5l The Derrida-Appert 
ratio {^t)c{^t)c/{^t)c froni [20] is also plotted in hgure O It has a non-zero hnite limit 
both when f —)■ 0 and f —)■ cxo, unlike the individual cumulants. 

4-3. Large deviations of the current in the long time limit 

At large time t, the generating function fHUH is essentially equal to the contribution of 
the stationary state corresponding to four empty sets, which implies 

(e*«t) ~ etuis) ^ ^43^ 

with 


fsfls) = xo(27rco(s)) 


(44) 
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Figure 6. Derrida-Appert ratio {^t)c{^t)c/{^t)c plotted as a function of the rescaled 
time t for an evolution conditioned on flat initial and final states, along with large and 
small t asymptotics from dUD and O- 


and co(s) solution of Xo(27rco(s)) = s. The function fst is the stationary state cumulant 
generating function characteristic of KPZ universality at large time. It corresponds 
to cumulants of the current scaling as {^^)c ~ t in the long time limit. It was hrst 
obtained for periodic TASEP in m by Derrida and Lebowitz, see also izni, and was 
subsequently derived for other models in the KPZ universality class: ASEP [ID], open 
TASEP on the transition line between low/high density phase and maximal current 
phase mi 021 Ha], discrete time ASEP with parallel update [Hji ^ directed polymer 
model [ID]) the asymmetric avalanche process [IDj. The stationary large deviation 
function was extended to the crossover between KPZ and equilibrium fluctuations in 
ASEP with weak asymmetry m SBl SHI EOj. Some results were also obtained for the 
average of current fluctuations in ASEP with several species of particles mm- 

At large time T, the total current Q for periodic TASEP is equal with probability 1 
to N{L — N)T/{L — 1) at leading order in T, see e.g. [5|. It implies —?■ 1 when t ^ oo 
for the current fluctuations. The Legendre transform gsflu) = maXs(Ms — /st(s)) of flt 
describes the probability of rare events when ~ tu with u 7 ^ 1 at large t. It is known 
as the large deviation function of the current for the stationary state, and verihes 

= tu) ~ . (45) 

In the notations of [20], one has fst{s) = G(s\/^)/\/^, gst{u) = —H{u)/y/Tr, 
and ^fst behaves for large argument as gst{u) — 2 \/ 3 M^/^/( 57 r) when u —)■ 00 and 
gsflu) ~ 4|Mp/^/3 when u —?■ — 00 . 

The functions fst and gst do not contain any information about the initial and 
the hnal state of the evolution. Some information about the evolution can however be 
found in the hrst order correction in t. In the hat —?■ hat case, one obtains from fldDj) 
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^ ^ ^ ^ , 27rco(s)-log\/^-log\/rTe2^^-logXo(27rco(s)) 

Jst [s) — Jst[s) H- - - 

up to exponentially small corrections in t. For the first cumulants of the current, it leads 
to 


(6) ^t+ 

2 \ y/2 VsJ 


(47) 




2^2 , 5 28 \ 3/2 




+ 


, . /15 9 24 \ 3/2 7319 20\/2 93 , ^ 96 \ 


71 


3 ' \/3 ' ^2 ^57 

These expressions are plotted in hgure |5] along with the exact hnite time values of the 
cumulants. 


4 . 4 . Large deviations of the current in the short time limit 

At short time, the hrst cumulants of ft scale as {ff)c ~ in the flat —?■ flat case, as 

seen in hgure O It corresponds for the generating function to the behaviour 

which is related by Legendre transform /o^^(s) = max„(sn — gl~^flu)) to the large 
deviations 

F(^/ = _ (49) 

This kind of short-time large deviations were already observed in [53] as a consequence 
of the fact that spatial correlations scale as for small t. They can be understood by 
breaking up the full system into around L/T'^^^ ~ t~‘^C almost stationary subsystems 
and by using stationary-like large deviations for each subsystem of size [53] . 

The function gl~^^ can be computed explicitly from a saddle point analysis of the 
exact formula (14^ for the probability at time f of see Appendix C One hnds 

/I /3 t-2/3S(«) 

t^^^Pflt^^^u) ~-=- (50) 


V^Zt 


with 


[u = 


dta w 


log(l 




E 

k=l 


Ai'(—2 ^/^/c^/^m) 

22 / 3/^573 


(51) 


in the range u G [0, (37r/2)2/3]. Outside of this interval, one has to add 3u/2 to 
the expressions for S in flSTD when -u < 0, and to add a.c.(u) dehned in (10.14^ for 
u > (37r/2)2/^. These additional terms make the function S analytic around the whole 
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Figure 7. Graph of the function S (EH), related to the large deviation function of the 
current at short time (1551) for an evolution conditioned on flat initial and final 
states. The asymptotics (1551) and (1551) are plotted with dashed lines. The dotted lines 
indicate the maximum of S. 


real axis. The function S is plotted in figure [3 The path of integration in required 
to avoid the branch cuts due to the logarithm plotted in hgure [8], goes to inhnity in 
directions specihed by angles ±6, ir/G < 6 < 7r/2. The symbol Ai denotes the Airy 
function. 

Calling n* 1.29146805131163785850008244580 the location of the maximum of 
S, the normalization condition = 1 implies the small t asymptotics 

for the ratio between the probability to observe the system in the flat conhguration iF 
at time t starting from F and the stationary probability of F. One has 


^2/3gt-2/3E:(u*) 


(52) 


This was checked numerically by truncating the sum over all eigenstates, see hgure O 
The short time large deviation function from fl49|) is then equal to 


slr'M = s(u,) - 3(u) , 


(53) 


with H(n*) = maXugRH(n) 0.360699035681939348898709742128. The functions gQ~^^ 
fl53|l and gst, Legendre transform of fl4T|l (see [20] for technical details about the required 
analytic continuation), are plotted in hgure [Till 

The location u* of the maximum of S is the deterministic limit of the random 
variable in the limit t ^ 0, and {ft) — Higher cumulants can be computed 
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Figure 8. Branch cuts of the integrand in the expression (EH for the quantity ^(m) 
for various values of the parameter u. The solid red curves represent the branch cuts 
of log (l — as a function of w. The dotted grey curves are the locus of the 

w such that Re(^ + 2'^Puw) = 0. The dashed black curves correspond to a suitable 
contour of integration oriented upward for EH- The different graphs correspond from 
left to right to u equal to —0.1, 0 (top row), 0.1, 3 (bottom row). 


numerically from fl53|) as well, and the leading correction obtained from ([50]), (|52|). One 
finds 


(ft/t) ~ 1.2914680513116378585 r^/^- 0.15751847779536355747 

^ 1.0126383442452675715 + 0.12793113905830599461 (54) 

^ -0.31901850108546379645^-2/^ - 0.19190248943176059601 
((^t/f)^), ~ 0.35959857102590971120 f-2/3 + 0.41299282842296155443 . 
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Figure 9. Normalization factor equal to the probability to observe the system 
in a flat configuration P for an evolution starting in P, divided by the stationary 
probability l/ff of P. The normalization is evaluated by summing over the 6639349 
eigenstates corresponding to sets A = Af = A~ and A = Aq = A+ with 
+ SogT® — quantity is plotted as a function 

of the rescaled time t (solid line). The dotted line is the limit t —>■ 0, equal to 
a; 1.006299331334999730981973520718. 



Figure 10. Graphs of shifted large deviation functions of the current as a function of 
their argument u. The solid red, lower curve represents the large deviation function 
at short time + u) for an evolution conditioned on flat initial and final 

states, defined in dsa). The solid blue, upper curve represents the stationary large 
deviation function gst(l + u) plotted from its Legendre transform (Hil) . The shifts of 
M* Ri 1.291468 and 1 are such that the minimum of both curves is located at 0. The 
black dots correspond to finite time evaluations at t = 2 (upper dots), t = 0.5 and 
t = 0.2 (lower dots) of —pP \og{pPp^{pP{u + u*))) from (1551) . after discretizing the 
integral over s and truncating the infinite sum over eigenstates in (I40|) . 
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These asymptotics are plotted in figure [5] along with exact numerical values of the 
cumulants. 

The asymptotics of the function S when its argument becomes large can be 
calculated explicitly. At large |m|, the expression flHT]) of S(m) is negligible compared to 
the extra terms 3u/2 and a.c.(u) required for the analytic continuation. When u ^ —oo, 
the asymptotics of the Airy function gives 


[u) ~ 


3^ (_^)i/4e-4(-.)^/V3 

2 2 \/^ 


(55) 


When M —)■ cx), the sum a.c.(u) flC.ldp becomes an integral, that can be computed 
explicitly. Adding the first Euler-Maclaurin correction, one finds 


[u) ~ 


U—^OO 


57r 4 


(56) 


We observe that the short and long time large deviations 5 fQ^^(M) and gsflu) have the 
same asymptotics ~ with the same coefficient in front when n —)■ cxo, see also figure 
m Similarities between short and long time large deviations were already observed 
from simulations in [53]. For u —oo on the other hand, gl~^flu) grows as |m|, slower 
than gst{u) which grows as 


4-5. Comparison with an evolution not conditioned on the final state: simulations 

Exact results for an evolution conditioned only on the initial state are still out of reach 
since they would require the large L asymptotics of the sum over all configurations 
for the first eigenstates, which is not known yet. It is however possible to 
study current fluctuations from simulations when the evolution is not conditioned on 
the final state. The first cumulants of the current are studied from simulations of 
periodic TASEP with N = 1000 particles on L = 2000 sites and fiat initial state. They 
are plotted along with the fiat fiat exact results in figure HH We use in this section 
the superscript f for the free evolution with fiat initial state, and the superscript 
f —>■ f for conditioning on fiat initial and final states. 

One finds the same scalings for the cumulants at short time in both cases, but 
with different rescaled cumulants. This is not surprising as one expects to have several 
universality classes at short time, similarly to what happens for KPZ universality on 
the infinite line [23]. On the other hand, in the long time limit corresponding to 
the stationary state, one finds the same cumulants (('Ct^O^)c — which are 

proportional to t for large t. A better correspondence at finite t is however obtained by 
decomposing the current fluctuations in the fiat fiat case as = ^o^t /2 + (t/ 2 ^t 
where ^o^t /2 and ^t/ 2 ^t represent the first and last half of the evolution. For large t, 
^o^t /2 and ^t/ 2 ^t are independent and have the same statistics as The additivity 
of cumulants for independent variables implies (('Ct^O^)c — 2((.^j(i^)^)c- Equivalently, 
— 2^“^ holds at large t, with rather good agreement 

for moderately large f, and very good agreement for the variance /c = 2 at small t too, 
see figure HH 
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Figure 11. Rescaled first cumulants of the current fluctuations plotted as a 
function of t. The solid lines correspond to for an evolution starting 

on the flat configuration and not conditioned on the final state, and are the result 
of simulations with N = 1000 particles on L = 2000 sites averaged over 326656 
independent realizations. The dashed lines correspond to for 

an evolution conditioned on flat initial and final states, computed from the generating 
function (HOI) . From top to bottom on the right side of the graph, they represent the 
rescaled average (blue), variance (red), fourth cumulant (orange) and third cumulant 
(cyan). The lower horizontal dotted line is times the mean value of GOE Tracy- 

Widom distribution. The upper horizontal dotted line is times the variance of 

GOE Tracy-Widom distribution. 


Apart from the total current Q, another interesting quantity is the (local, time- 
integrated) current Qi between sites i and i + 1 (at half-hlling p = 1/2 only because of 
the necessity to consider a moving reference frame with velocity 1 — 2 p in order to see 
fluctuations characteristic of KPZ universality). When the evolution is conditioned on 
both the initial and hnal conhguration, Q and Qi are closely related since particles 
can not overtake each other (this can also be understood by a simple similarity 
transformation of the deformed Markov matrix M{'-y) with a diagonal change of basis 
in the conhguration basis). In particular, if the initial and the hnal conhgurations are 
identical, one has Q = LQi at the hnal time T. This is not the case any more for an 
evolution conditioned only on the initial conhguration. 

On the inhnite line Z, the statistics of Qi have been investigated in much detail 
[21] for an evolution not conditioned on the hnal state. For hat initial condition, the 
probability density of with = 2(Qj^ — JT)/\/X and ffT7|l is given 

[511155] by the (derivative of the) GOE Tracy-Widom distribution from random matrix 
theory. This result presumably also holds at any rescaled time t = T/L on the Euler 
time scale T L, and in the limit t —?■ 0 on the KPZ time scale. If the initial condition 
has particles only at odd sites, we observe rather strong hnite size corrections for the 
mean value at small t of /t^^^ with i even, which disappear for i odd, see hgnre [T2l 

The precise relation between the GOE Tracy-Widom distribution for the 
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Figure 12. Mean value (left) and variance (right) of current fluctuations for periodic 
TASEP with flat initial condition and no conditioning on the final state, plotted in 
terms of the rescaled time t CZl). The solid lines, corresponding to the mean value and 
variance of the fluctuations of the current Qi between sites i and i + 1, 

are the results of simulations with N = 1000 particles on L = 2000 sites, averaged over 
134144 realizations for an initial condition with particles only on odd sites. On the left, 
the upper solid line, which shows large finite size corrections for small t, corresponds 
to * = 0 while the lower solid line, which does not show any noticeable finite size 
effect, corresponds to i = 1; on the right, the two solid lines are almost superposed. 
The dashed lines correspond to the fluctuations of the total current Q. On the 
left, the mean value —is plotted, and on the right the rescaled variance 
2i/3^-4/3(^2)^^ The dotted lines correspond to the mean value (left) and the variance 
(right) of GOE Tracy-Widom distribution. 


fluctuations of and the large deviations of the total current is currently not known. 
However, since ■, their mean values must be equal. More precisely, for 

an initial condition with particles only at odd sites, the mean value of ft is equal to 
for a finite system. The hnite-size corrections to for small t are 
responsible for the not so good convergence to the mean value of GOE Tracy-Widom 
in figure [HI see also hgure [121 We have no explanation however for the numerical 
coincidences for the variance at small time 

1.607781034581 ^ VarcoE = 

^lim2V3t-4/3((^f-)2^^ ( 57 ) 

liin22/3r^/3((^f^f)2)^ ~ 1.6074631729182734577 , 

with VarcoE the variance of GOE Tracy-Widom distribution, see f[5T|) . figure [H] and 
figure [HI 

5. Conclusions 

Gurrent fluctuations for periodic TASEP on the relaxation scale are studied in this paper 
using large system size asymptotics of eigenvalues and eigenvectors of the generator of 
the evolution. For technical reasons, our results are restricted to evolutions conditioned 
on simple initial and final states. An exact formula for the generating function of current 
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fluctuations is obtained as a sum over eigenstates. In the special case of flat initial and 
final configurations, it leads to a simple expression fl3T]) for the large deviations of the 
current at the early stages of the relaxation, written in terms of the Airy function. 

Extending these results to more general initial and hnal states would be interesting 
in order to fully describe the process on the relaxation scale. Removing the conditioning 
over the hnal state would also allow to understand better the relation with the Tracy- 
Widom distributions that describe current huctuations on the inhnite line, by hnding 
short time large deviation functions for the total current taking their minimum at the 
mean value of the corresponding Tracy-Widom distribution. 

The results obtained in this paper should presumably extend to all models in one¬ 
dimensional KPZ universality. It would be interesting to recover them directly from 
stochastic Burgers’ / KPZ equation with periodic boundary conditions, using the replica 
method with precise asymptotics for the attractive (5-Bose gas in hnite volume. 
Acknowledgements: It is a pleasure to thank C. Bahadoran, M. Bauer, K. Mallick 
and H. Spohn for useful discussions. 

Appendix A. Total integrated current for periodic Burgers’ equation 

In this appendix, we check the expression ffT5|) for the hrst correction at large times 
to the total current in the inviscid Burgers’ equation (|1]) by working out a few simple 
examples with piecewise linear initial condition po, for which the solution of the partial 
differential equation can be computed easily: the linear parts evolve in time as ^ 
and a discontinuity at xq with densities respectively p_ and on the left and on the 
right leads if p_ < to a shock that moves with a velocity 1 — p+ (r) — p_ (r) imposed 
by conservation of density, and if p_ > p+ to a rarefaction fan with density prohle 
p{x, t) = \ — in fhe interval xq -l- (1 — 2p_)r < x < xq -l- (1 — 2p+)r until one of 
the extremities is absorbed by a shock. 

Appendix A.l. Step initial condition 

The step initial condition is dehned as the periodic prohle with periodicity 1 equal 
for X e [0,1) to po(x) = p+l{ 3 ;e[ 0 ,a)} -f P-l{a;e[a,i)}, with 0 < p_ < p+ < 1 and 
0 < a < 1. The corresponding average density is p = ap+ -|- (1 — a)p_. We 
assume a < ^, the case a > ^ being accessible by changing po to 1 — po. At the 
beginning of the evolution, the shock initially located at zq = 0 moves at the constant 
speed 1 — p+ — p_, while a rarefaction fan opens at position a, with its left side 
moving at speed 1 — 2p+ and its right side at speed 1 — 2p_, see hgure lAll This 
continues until time Ti = a/(p_|_ — p_) when the interval of density p+ disappears as 
the left of the fan merges with the shock. Then, the position of the shock moves as 
z{t) = a — 2^/ay/p+ — P-y/r -|- (1 — 2p_)r until time T 2 = (4a(p+ — P-))“^ when the 
interval of density p_ disappears as the right of the fan merges with the shock. Finally, 
after r 2 , the shock moves as z{t) = (1 — 2 p)r + a — The corresponding density prohle 
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Figure Al. Evolution of the density profile from Burgers’ equation (U) with periodic 
boundary condition for unit step (left) and linear increasing (right) initial condition at 
average density p = 1/3. The horizontal axis correspond to space and the vertical axis 
to time. The solid, red curve represents the shock. Dashed, blue lines divide regions 
with a different behaviour for the density profile. Dotted horizontal line indicate the 
times at which edges of rarefaction fans merge with the shock. 


is then p{x,t) = | — for z{t) < x < z{t) + 1. Equivalently, in a reference frame 
moving at velocity 1 — 2p, one has p{x + (1 — 2p)r, r) =p — for a — ^ < x < a + |. 
After some calculations, one hnds for the total integrated current up to time r > T 2 


QtIpo] = p{^- p)r 


{p- P-){p+- p) 


+ 


(A.l) 


2(p+ - P-) ' 48r 

The constant term agrees with the general expression flT^ for TZ[po] with k = a — 1 
modulo 1. 


Appendix A.2. Linear decreasing profile 

We consider the initial prohle po{x) = p — a{x — ^) for 0 < x < 1 and 0 < a < 
2 min(p, 1 — p), with a shock initially located at zq = 0. At any time r > 0, the suitable 
solution of dH) is equal in the moving frame to p{x + (1 — 2 p)r, r) = p — 

0 < x < 1. This implies k = ^ and a position of the shock z(t) = (1 — 2p)r. The 
calculation of the current from the previous expression for p(x, r) gives 


Qr[po] =p(l - p)t- 




12 + 24Q!r 


(A.2) 


_x « I 


I 


+ 


96ar^ l92a‘^T^ 


The constant term matches with flTHli . 
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Appendix A.3. Linear increasing profile 

We consider the initial profile po{x) = p + Q!(x — ^) for 0 < x < 1 and 0 < a < 
2niin(p, 1 — p). In the beginning of the evolntion, a rarefaction fan opens at position 
0, see hgnre lATl Until time ti = the density prohle is p{x,t) = | — ^ 
for X in the interval [(1 — 2p + a)T, (1 — 2p — a)r + 1] and p{x, r) = | ^ for x in the 

interval [(1 — 2p — a)r, (1 — 2p + a)r]. Then, at time ti, the linearly increasing portion 
vanishes and a shock forms leading to the density prohle p(x + (1 — 2 p)r, r) = P — ^ 
for — ^ < X < ^. A space-time integration leads for r > ti to 

+ (A,3) 

which can be recovered directly from ffTT|) . flT5|) with k = 0. This expression is however 
exact for r > Ti, nnlike flTT|) which is the beginning of a large r asymptotics. 


Appendix B. Bethe ansatz for TASEP 

In this appendix, we summarize some known results about Bethe ansatz for the 
hrst eigenstates of TASEP, in particular large L, N asymptotics of normalization of 
eigenvectors and components of the eigenvectors corresponding to unit step density 
prohle. We also derive the asymptotics fl33l) for the components of the eigenvectors 
corresponding to a hat density prohle. 


Appendix B.l. Bethe eguations and their solution 


From Bethe ansatz, each eigenstate of periodic TASEP for a hnite system of length L 
with N particles is completely characterized by N complex numbers yj, j = 1,... ,N, 
the Bethe roots, that satisfy a set of N polynomial equations called the Bethe eguations: 

N 

e^'>(l-p,)^ = (-1)^-1 J]-- (B-1) 

Multiplying both sides of the Bethe equations by and taking the power 1/L 

gives |56l 130] 

2i7rfc.,- 

g{yj) = e—, (B.2) 


where the kjA, distinct modulo L, are integers (half-integers) if N is odd (even). The 
function g is dehned by 


9{y) = 


and b is solution of 


1 - 2 / 


N 




i=i 


(B.3) 


(B.4) 
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The branch cut of g due to the non-integer power is taken as which leads to the 
branch cuts cxo) for the inverse function g~^. 

The equation fIB.dp can be solved numerically for b with high accuracy using 
Newton’s method. At each step, the yfls are computed by inverting g in flB.2p using 
again Newton’s method. It is possible to obtain very accurate expressions for the Bethe 
roots i/j, with several hundred significant digits. Such accurate values are needed in order 
to fully exploit the power of Richardson extrapolation for obtaining precise asymptotics 
of various quantities from a few hnite size values, see table [U for an example. 

The question of the completeness of Bethe ansatz for hnite systems has not been 
fully solved yet, see however [571 [58l |59l [56] . Nevertheless, one observes that the number 
of possible choices for the kfls with 0 </ci<...</cAr<Lis equal to the number of 
conhgurations Q = (^). Comparison with exact diagonalization for systems up to size 
L = 18 seems to indicate that each such choice of the kfls corresponds to an eigenstate 
of the Markov matrix with flB.dp having a unique hnite solution, if Re 7 is large enough 
so that the solution b of flB.dp is such that the circle of center 0 and radius e“^ does 
not cross the branch cuts of g~^. Numerics indicate that Rey > 0 is sufficient for all 
the eigenstates. The stationary state kj=j — {N + l)/2 with 7 close to 0 is special 
since the circle crosses the branch cuts of g~^, which is however not a problem since the 
points stay on the same side of the branch cuts in that case; this special case 

leads to the singular solution co(s) —)■ —cxo when s —)■ 0 of fl26p in the large L limit. 

An alternative approach was used in [TSl Hn] to characterize the solutions of the 
Bethe equations fIB.ip at 7 = 0, by rewriting gigjj)^ = e“^'^ as T’((l +yj)/{l —yfl) with 
P a polynomial of degree L with coefficients depending on b. Using a particular labelling 
of the L roots of P, the eigenstates were then identihed as choices of N distinct roots 
of P among L. As noted in [13], however, this identihcation fails for some eigenstates 
of large enough systems. This is due to the fact that with this specihc labelling for the 
roots of P, changing the imaginary part of b can induce a cyclic relabelling of the roots 
of P. 

Appendix B.2. Eigenvalues and eigenvectors 

The eigenvalue of M ( 7 ) corresponding to a given solution of the Bethe equations fIB.ip 
is given by 

N 

= E ■ (B-5) 

i=i 

The corresponding eigenvalue for the translation operator U, U\xi,... ,xisf) = |1 + 
xi,..., 1 -|- xn), is equal to 

N N 

^2inpr/L ^ _ y.^ ^ ^2ink,/L _ 

i=i i=i 

The components of the left and right eigenvectors of M( 7 ) with particles at positions 
Xj, j = with 1 < xi < ... < xat < L are given by the Bethe ansatz as 
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linear combinations of all N\ permntations of N plane waves with pseudo-momenta 
7 -|- log(l — yk). For TASEP, the sum over permutations reduces to a determinant. One 
has 




JV+l 


i=i 

N 


I X 

(^.|x) = i^e^AT-^/^ ( ) det fe(l - 




i=i 


j,k=l,...,N 


(B.7) 
. (B.8) 


These Bethe eigenstates are not normalized. The factors in front of the determinants 
are chosen in prevision for the thermodynamic limit. 

For any conhguration x, the reversed conhguration x is dehned hy Xj = L + 1 — 
xjY+i_j. One has the symmetry relation 


(^|x) = {Stlfj) , (B.9) 

which is a consequence of the fact that transposing the evolution operator of TASEP is 
equivalent to reversing space. 


Appendix B. 3. Normalization of Bethe eigenstates 


The norm of Bethe eigenstates is in general given by the Gaudin determinant 
It reduces for TASEP to the explicit expression 




N . 

{TT 


V j=i 


+ Vj ^ 

P ^ 


(B.IO) 


whose asymptotics can be obtained using the Euler-Maclaurin formula. With O = 
the total number of conhgurations, one has [2B] 

n _ (27r)-F^e^""- (na6A+ + iQ)(nagA- - la) 

\/rTe^x"(27rc^) (naeA+ \/Cr - ia)(naeA- Vcr + ia) ’ 
where Xr is dehned in fl25|) and Cr is the solution of fl26|) . This leads to fl29|) after changing 
the normalization of the eigenvectors as (0^1 = A“^('0r| and |0r) = A“^|'0r) with 






- “)) 


(1 +e^""’') 


27rCr\l/4 


X 


n 


a&A+ (g - ia) ( riaeA- (g + ia) 


( naeA+(cr + ia)^/^)( riaeviy (cr - 


(B.12) 


Appendix Bfl. Flat configuration (1/p integer) 

We consider the hat conhguration IFx with particles at positions xj = X -|- (j — 1) /p, j = 
1,... ,N and p~^ integer. The determinants in flB.Tp and fIB.Sp are then Vandermonde 
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determinants: 

N{N-1) TnrprX {N-l)L-i 


N 




x\ 


= -1 


e L e 


Hi'. 


N-1 

2 


j = l 


(B.13) 


l<j<l<N 

where the function g is defined by flB.3p . This expression can be simplified further by 
noting that the Bethe equations precisely give an explicit expression flB.2p for g{yj) in 
terms of the (half-)integers kj. Using flB.4D to simplify the single product of the yfls, 
one has 

. (b.14) 

l<j<l<N 

From the symmetry relation flB.Qp . the left eigenstate is given by (fjliFx) = {J^L+i-x\fl’)- 
We observe that these expressions are non-zero if and only if the kfls are all distinct 
modulo N. For the first eigenstates, described in figure [21 it is equivalent to the 
constraints Aq = A~ and Aq = A~^, which imply that all the sets have the same cardinal 
nir = rnf = m~ and that the total momentum pr = 0. Splitting the contributions to 
the double product coming from the Fermi sea k^ = j — {N -|- l)/2 and from the sets 
At, one hnds 


a)+(E_.- o) 


JJ" ^g2i7rfc^/Af _ g2i7rfcj/V^ _ + JJ" ^g2i7rfcO/Af _ g2i7rfcP/V^ ^ 

The double product reduces to a simple product by factoring out and making 

the change of variables j j + i. The remaining simple product can be computed using 
the symmetry j N — j. One hnds 

J-J ^g2i^fc«/JV _ g2i^fc°/Af) = •^4 ^^jv/2 _ 

It hnally gives the asymptotics with independent of X, 

given by 


which leads to fl55p after the change of normalization above flB.12p . 


(B.17) 


Appendix B. 5. Step configuration 


We consider the step conhguration Sx with particles at positions Xj = X + j — 1. The 
corresponding component of the eigenvectors are 


(5 


JV(iV-l) 2iTs-prX N{N-l)-t 


N 


N N 




= -1 


^ "''diis'j" )n n (»-»). 


(B.18) 


j=l j=l k=j+l 

and {iplSx) = {SL-N+ 2 -x\fl’)- The large L asymptotics of this expression was studied in 
[28] using two-dimensional Euler-Maclaurin formula with various logarithmic and square 
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root singularities at the borders of the summation range. It has an expansion in powers 


of l/y/Tj instead of 1/L for the flat case. One has 



(B.19) 


with 



(B.20) 


(27r)”^’- 


2-KCrVr/i ( nagA+ [Cr - iaV^^){Ua&A-iCr + 






The combinatorial factors u are dehned in fl5B]) . This leads to fl5SD after the change of 


normalization above flB.12p that cancels the third line and some factors in the hrst line 


of (IR^ . 


Appendix C. Saddle point analysis of the flat —flat case at short time 

In this appendix, we derive the expression fl3T]) . for the large deviation function 
of the current at short time for an evolution conditioned on flat initial and flnal 
conflgurations. 

We start from the exact formula flT2l) for the probability density of ft at arbitrary 
rescaled time t, and consider instead the probability density of Making the 

change of variables c = t~‘^^^d and using 



Rex > 0 
Re X < 0 


(C.l) 


the small t limit of the integrand gives 



(C.2) 


with 



(C.3) 


The branch cuts of the square roots are chosen equal to so that \/n ± id is analytic 
in both half-planes Red > 0 and Red < 0. Because of the branch cuts i[|,cx)) and 
— i[|, cxd) in the variable c, the two half-planes become independent when f —)■ 0. 
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It is convenient to make the change of variables w = 2^/^\/i\/4^\/n + id and 
w = — id respectively in the hrst and the second integral. This leads to 


h(^U,djZ) l{Re(i>0} 

w dtc 


-u -h 71 d 


2071 


Stt 


+ 


'r+ 


25/3i 


ITT 


log 


1 + 


+ l{Re(i<0}27r(i 

w dtc 


(C.4) 




The contour r_|_ is oriented from a+ = 2^/^\/i\/4^\/id to \/ioo, while the contour r_ is 
oriented from \/^cxo to a_ = 

In Appendix C.l and Appendix C.2 we show that, for d G with some specihc 
unbounded domain of the complex plane, the relation ddh{u, d,z) = i log(—holds, and 
the quantity h{u, d, z) is thus affine in d. The saddle point equation for the variable d 
then gives z = —1. Expanding at second order in z, one has 


h{u, d, z) ~ S(m) — i{d — d^,(u)){z + 1) + {A{u) — i(d — dflu))) 


{z + lfl 


+ 0{z + Ifl ,{C.5) 


with H(m) = h{u, d, —1) independent of d, and some functions (i*(n) and A{u). In order 
to perform the Gaussian integration around the saddle point, we write z = —1 + aiy 
and d = dflu) + uix, a = ±1, z/ = ±1 (it does not matter any more at this point 
whether dflu) still belongs to the domain since, h{u, d, z) being affine in d, it can be 
continued analytically to all d G C before moving the contour for d through the saddle 
point dflu)). The integration over y gives 


Zt 


I t-^73x^ \ 

dx ^M- 2iAiu}+..) ) 

27r \/ A{u) + nx 


(C.6) 


The saddle point for x is x = 0 and one can then replace A{u) + vx by A{u) in the 
integrand. After Gaussian integration in m, one hnds (15011 with 5 (m) = h{u, d, —1). The 
signs have to be chosen as az/ = — 1. 

For technical reasons related to the complicated branch cut structure of the 
integrands in (IG.4I1 as a function of w, see hgure [3 the choice of a suitable domain 
Vu for d depends on the sign of u. When m > 0, we will need to consider a domain 
included in the half-plane Red > 0, while for u < 0, the domain will be included 
in the half-plane Red < 0. As shown on hgure 01 it is always possible to deform the 
contours Cj.(is), s G M in the hnite time expression fl55]) . (00]) to a unique contour 
with Rec > 0, which makes the saddle point analysis straightforward when u > 0. 
Deforming the contours Cr(is), s G M to a contour with Rec < 0 is on the other hand 
not possible. One can however always make some portion of the contours pass through 
the hole i[—between the branch cuts in hgure 01 The hole then closes after making 
the change of variable from c to d and taking the limit t —?■ 0, and the integral over d in 
(IG.2D can be decomposed as an integral with Red > 0 plus an integral with Red < 0. 
Since only the part with Red < 0 seems to possess a proper saddle point, we assume 
that the contribution of the integral with Re d > 0 is negligible in the small t limit when 
M < 0. This seems justihed by the fact that the expression obtained in the end for the 
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Figure Cl. Domain Vu defined in (1C.91) for m > 0 and in (1C.171) for u < 0, plotted 
from left to right for u = —3, u = —2, u = —1, u = 1, u = 3 and it = 5. 


large deviation function is analytic in u, and by comparison with numerical evaluations 
of the probability density at small times, see hgure [Till 


Appendix C.l. Case m > 0 

We consider here that 2 ; and d, Red > 0 are such that it is possible to choose a 
determination of the logarithm so that its branch cut is never crossed in fiC.djl . Taking 
the derivative with respect to d of fiC.djl and using log(l + q~^) = log(l + q) — logg for 
q ^ , one finds for generic values of d 

ddh{u, d, z) = + TT + i \og{-z e-+f . (C.7) 


Writing z = r e‘^, r > 0, 6> G M, it leads to 


ddh{u, d,z) = i log(- 2 ) + 271 


27r + d + Im(|^(47rd)^/^ — 2m\/47rd) 


27r 


(C.8) 


with [xj the largest integer lower than x. We define a domain C C, equal to (the 
unbounded connected component of) 

27r + Im(y (dvrd)^/^ — 2m\/47rd) 


Vu = {d, Re d > 0, 


27r 


= 0 }. 


(C.9) 


The domain always contains a suitable path between and e see figure 

EH 

For d E Vu and z in a neighbouring of —1, one has ddh{u,d,z) = ilog(—z), hence 
the function S(m) = h{u, d, —1) is independent of d. Its expression obtained from flC.411 
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Figure C2. Contours r±, A+, A_ and —A+ (blue, green, yellow, cyan, dashed lines) 
for M = 5 and d = 1.1 G plotted along with the branch cuts (red, dotted lines) 
and poles (red dots; after integration by parts) in the variable w of the integrand in 
(IC.13I) . All contours are oriented upward except — A+ that is oriented downward. 


is equal to the sum of two integrals with the same integrand, and contours of integration 
r+ and r_. Since Re d > 0, the finite ends a+ and a_ of the contours verify a++Q!_ = 0. 

The expression of S(n) from flC.djl can be simplihed by closing the contour between 
a_ and a+ on the path A+UA_ represented in figure for generic n > 0 and d G 
(which implies that a± is purely imaginary, and the contours r± stay in the half-plane 
with positive real part; if d ^ M"*", either Ima+ < 0 or Ima_ <0 and the branch cuts 
with Retc < 0 have then to be chosen in such a way that they do not intersect the 
contours r±). After an integration by parts to replace branch cuts with poles, using 
again log(l -b q~^) = log(l + q) — logg (for generic q) and flC.9p . one finds 


'A_UA+ 


w dtp 

2^ 


ITT 


log 1 


gif+22/3^ 


(dtrd)®/^ (dvrd)^/^ 


Gtt 


u ■ 


27r 


-b vrd (C.IO) 


tp^ dtp tp^ -b 2 ^/^m 


ITT 


1 _ 


dA_uA+ 

All branch points have become poles except tp = 0 for which the residue vanishes. 


Replacing w by — tP in the integral on the contour A_|_ and using (1 — e’s 


■+ 2 '^/^uw\-l _ 


l-(l-e 


-22/3uu;^-l 


) the integral in the right hand side of fIC.lOp becomes an integral 


over the closed path A_ U (—A+), see figure (221 We distinguish two kinds of poles inside 
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this contour: ’’inner poles” wfl, k = 1,... ,n, located closest from 0, that the contour 
encircles clockwise, and ’’outer poles” k = 1,... ,n, located farthest from 0, that 

the contour encircles counter-clockwise. The integer n is equal to the number of strictly 
positive k such that ^ -|- 2‘^/^uWk = —2mk has 3 real roots iwk- This is equivalent to 
A = 16n^/3 — 3(27r/c)^ > 0, where A is the discriminant of the third degree equation. 
One has thus n = • The residues of inner and outer poles are equal to 


Res 

w=Wf^ 2 ! 


w 

Res 


_|_ 2‘^/^u 

_ g-3^-22/3n, 

yj2 _|_ 2'^/^u 


= —2Msin^ 


= —2m sin 


arcsin ^ 


-TT arcsin 


^ _ ^-^-2^1-^UW \ 3 

One hnally hnds that the integral on the path from a_ to is equal to 


'■“+ wdw 


25/3i 


ITT 


/ A 

log (1 — e 3 


■+22/3mw,\ _ 


(47rd)5/2 

ISvr 


(C.ll) 

( 0 . 12 ) 


■M^^^5^ + 7rd-a.c.(M) ,(C.13) 
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with 


1^1 
L Stt -I 


a.c.(M) = 2u 


k=l 


sm 


arcsin 


Snk 

2«3/2 


sm 


-TT P arcsin 


. (0.14) 


This leads to the integral expression in flSip for the function S(m). We observe that the 
additional term a.c.(M) is precisely the one required to make S analytic near the positive 
real axis. 

The function S(m) can also be expressed in terms of the Airy function by expanding 
the logarithm and using 


Ai(g) 


r‘icxD+e 


' —icx)+e 


dtc 

2i7r 


e 



(0.15) 


with e > 0. This leads to the expression with the inhnite sum of S in fl5T]) . 


Appendix C.2. Case m < 0 


We consider here d in the half-plane Red < 0. It implies that the finite extremities a± 
of the contours r± are equal. The expression fl0.4p of h{u,d, —1) can then be written 
in terms of a single integral with a connected path of integration T = r_|_ U r_ between 
\/^cxo and \/icxo: 


h{u, d, —1) 


We introduce the domain 


V 


U 



3^2 


16\/3u — dvr Red 


< Re d < 


3m 'i 
dvr/ 


5 


(0.16) 


(0.17) 


which is a vertical strip in the complex plane, see figure 1011 Using Re(M;^) = 
(ReMi)^ — 3(ReM;)(ImM;)^, Im(M;^) = —(Imtc)^ -|- 3(ReM;)^(ImM;) and (Rewfl = 
(ImM;)^ — 2^/^7rRed for m; G T, one can show that for d G the contour T crosses the 
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Figure C3. Contour F = r+ ur_ (blue, dashed line) for u = —1 and d = —0.3 € 
plotted along with the branch cuts (red, dotted lines) and poles (red dots; after 

integration by parts) in the variable w of the integrand in (IC.16|) . The contour is 

oriented upward. 


real branch cut (an only that branch cut) of the integrand in the expression fIC.lbp of 
h{u, d, —1), see hgure IC3I 

Integrating flC.1611 by parts, we observe that the contributions of both sides of the 
branch cut cancel the term 2nd. We obtain 

dtp + 2^/^u 


h{u, d, —1) = — 




(C.18) 


2^Ci7i p 

The point 2 ; = —1 is again a saddle point with respect to d, ddh{u,d,—l) = 0. 
Moving the contour to the left, one picks the term 3 m/ 2 from the residue of the 
pole wq = 2^/^\/3\/— M. Integrating by parts again, we recover 6D with a contour 
of integration that does not cross any branch cut. Again, we observe that the extra 
term 3 m/ 2 corresponds precisely to the analytic continuation to m < 0. 
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